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It  is  well  known  that  statistical  estimates  from  simulation,  e.g.  computation  of 
mean  and  standard  deviation,  often  involve  significant  uncertainty  due  to 
uncertainty  in  the  input  parameters.  The  uncertainty  in  the  output  results  can  be 
reduced  with  additional  data  on  the  inputs,  e.g.,  additional  experiments.  An 
optimization  methodology  is  proposed  and  implemented  to  determine  the  optimal 
additional  experiments  to  conduct  in  order  to  minimize  the  variance  in  the  output 
mean  and  standard  deviation  subject  to  a  cost  constraint.  The  number  of  additional 
experiments  to  add  for  each  random  variable  depends  upon  several  factors:  the 
number  of  initial  data  points,  the  importance  of  the  random  variable  in  the 
response,  the  range  of  the  random  variable,  and  the  cost  of  each  experiment.  The 
methodology  is  demonstrated  using  several  numerical  examples.  The  results  indicate 
that  the  particle  swarm  optimization  performs  well  and  the  solutions  obtained  are 
superior  to  ad  hoc  ones. 
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1.  Nomenclature 

=  Constants  in  the  response  function  for  variable  X. 

=  Resource  available  for  additional  experiments 
=  Cost  of  each  additional  experiment  for  variable  X. 

=  Number  of  additional  experiments  for  variable  X. 
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Number  of  initial  data  points  for  variable  X. 
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Optimization  iteration  number 

Lower  bound  at  optimization  variables 

Monte  Carlo  Sampling 

Number  of  total  data  points  for  variable  X 

Number  of  variables  in  the  response  function 

Personal  best  of  the  jth  pso  particle  until  the  current  iteration 

Individual  weight 
Social  weight 

Random  number  between  0  and  1  for  jth  pso  particle  and  first  term 
Random  number  for  jth  pso  particle  and  second  term 
Sample  standard  deviation  for  variable  X- 
Position  of  jth  pso  particle  at  k  iteration 

Student’s  t-distribution 
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Velocity  of  jth  pso  particle  at  k+1  iteration 

Inertia  weight  of  the  jth  pso  particle 

Values  of  X  for  Ei  data  points 

Values  of  X i  for  D  additional  experiments 

Statistically  independent  variable  in  the  response  function 

Sample  mean  for  variable  X. 

Response 

Chi-square  distribution 
Population  mean  for  variable  X. 

Response  mean 

Population  standard  deviation  for  variable  X- 
Standard  deviation  of  the  response  mean 

Standard  deviation  of  the  response  mean  based  on  the  original  data 
Optimum  solution  of  the  standard  deviation  of  the  response  mean 
Response  standard  deviation 
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1.  Introduction 


It  is  well  known  that  statistical  estimates  from  simulation  often  involve  significant  uncertainty  due 
to  uncertainty  in  the  input  parameters,  e.g.  the  probability-of-failure  or  the  mean  and  standard  deviation  of 
the  output.  The  input  probability  distributions  may  be  developed  from  limited  data,  e.g.,  a  few  experiments, 
thus  the  mean  and  standard  deviation  of  the  inputs  are  also  random  variables  dependent  upon  the  amount  of 
data  available.  The  resulting  uncertainty  in  the  simulation  outputs  can  be  very  large  and  should  be 
accounted  for  in  reporting  the  computational  results  and  subsequent  decision  making. 

Other  authors  have  recognized  this  issue  and  presented  strategies  for  determining  the  variation  in 
the  output  results  due  to  the  uncertainty  in  the  input  PDF  parameters.  While  the  concept  to  accomplish  this 
is  straightforward  by  nesting  the  single  output  probability  calculation  (inner  loop)  inside  another  loop 
which  accounts  for  the  variation  in  the  input  parameters  (outer  loop),  the  practical  accomplishment  of  such 
nested  probability  calculations  is  nontrivial.  Several  authors  have  developed  computational  strategies  to 
address  this  issue.  The  strategy  to  date  usually  involved  approximating  the  inner  loop  using  a  surrogate 
model  (response  surface),  then  calling  the  surrogate  model  from  the  outer  loop.  Thus,  the  analysis  is 
efficient  once  the  response  surface  has  been  constructed.  However,  the  quality  of  the  final  solution  will 
depend  upon  the  quality  of  the  surrogate  model. 

Cruse  and  Brown  present  a  Bayesian  analysis  method  to  compute  the  confidence  interval  of  the 
probability-of-failure  using  a  directed  acyclic  graph  approach  [6],  The  important  random  variables  were 
first  determined  using  sensitivity  analysis  (4  total),  then  these  variables  were  used  in  a  design  of 
experiments  study  to  construct  the  response  surface.  The  conditional  probability-of-failure  as  a  function  of 
the  random  variable  input  parameters  was  approximated  using  a  polynomial  response  surface  method  using 
design  of  experiments. 

Brown  and  Grandhi  develop  a  response  surface  of  the  probability-of-failure  with  respect  to  the 
input  PDF  moments  [2],  A  derivative-enhanced  Kriging  model  is  used  to  develop  the  response  surface.  The 
computational  complexity  issue  is  then  circumvented  by  replacing  the  inner  loop  with  the  response  surface. 
The  results  show  the  Kriging  model  performs  very  well  and  confidence  intervals  of  the  probability-of- 
failure  can  be  estimated  in  reasonable  time.  The  results  also  show  that  the  variation  in  the  input  parameters 
can  cause  dramatic  variations  in  the  probability-of-failure,  resulting  in  very  broad  95%  confidence  limits. 

Ferson  and  Ginzburg  emphasize  that  the  nested  approach  (called  second-order  Monte  Carlo  in 
their  paper)  has  inherent  difficulties:  (1)  parameterization  issues,  (2)  computational  complexity,  and  (3) 
interpretational  difficulty  [5],  Parameterization  requires  the  analyst  to  select  how  the  parameters  of  the 
inputs  vary  which  may  introduce  more  uncertainty.  However,  in  certain  circumstances  the  distribution  of 
the  input  parameters  is  well  known  theoretically.  For  example,  the  mean  follows  a  t  distribution  and  the 
standard  deviation  follows  a  chi-squared  distribution  [1],  Computational  complexity  is  a  continuing  issue 
that  requires  clever  strategies.  The  best  approach  to  date  is  to  develop  a  surrogate  model  of  some  form  for 
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the  inner  loop  that  defines  the  probabilistic  output  as  a  function  of  the  parameters  of  the  input  PDFs.  The 
interpretational  approach  recommended  here  is  to  compute  confidence  intervals  (typically  95%). 

While  several  authors  have  recognized  the  issue  of  uncertainty  in  the  input  parameter  PDFs  and 
the  computational  challenges  that  arise,  less  work  has  been  done  on  assessing  the  resulting  “action  items” 
needed  to  reduce  the  confidence  limits  on  the  probability  output(s).  In  this  paper,  an  optimization  method  is 
combined  with  a  simplified  inner  loop  scheme  to  determine  the  subsequent  experiments  needed  to  most 
efficiently  reduce  the  output  variance  subject  to  a  cost  (time,  money,  schedule)  constraint.  The  action  items 
considered  in  this  work  were  to  add  additional  data  (experiments)  to  better  characterize  the  mean  and 
standard  deviation  of  the  input  PDFs.  In  other  words,  what  additional  experiments  should  be  conducted  in 
order  to  improve  the  confidence  in  a  probabilistic  solution  given  cost  constraints. 


2.  Methodology 


A  schematic  of  the  computational  approach  is  shown  in  Figure  1.  The  inner  loop  is  shown  simply 
as  “Probabilistic  Analysis”.  The  outer  loop  indicates  that  a  Monte  Carlo  sampling  method  will  be  used  to 
generate  realizations  of  the  mean  and  standard  deviation  of  each  random  variable,  evaluate  the  mean  and 
standard  deviation  of  a  response,  then  repeat  to  develop  statistics  of  the  mean  and  standard  deviation  of  the 
output  response. 

From  statistical  theory,  the  quantity  ( X  -  fX)/{S I'ln)  for  an  independent  random  variable  has  a 
Student’s  t  distribution  with  n  -1  degrees  of  freedom,  where  X  denotes  the  sample  mean,  fl  the  true 
population  mean,  S  the  sample  standard  deviation,  and  n  the  number  of  data  points  [1],  As  a  result, 
realizations  of  the  population  mean  for  random  variable  X,  that  represent  the  expected  variation  in  the 
population  mean  as  a  function  of  the  number  of  data  points  are  straightforward  to  determine  as 


jU.=Xi  +  %^L 


(1) 


where  t  represents  a  realization  from  the  standard  t  distribution. 

Similarly,  the  quantity  (n  -  1)S2  /cr  follows  a  chi-square  distribution  with  n  -  1  degrees  of 
freedom.  A  realization  for  the  population  standard  deviation  can  be  determined  from 


Xn- 


(2) 


where  x  is  the  chi-squared  distribution  with  n  - 1  degrees  of  freedom. 

A.  Probabilistic  Analysis 
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The  determination  of  the  probability-of-failure  or  the  response  moments  is  often  computed  using 
sampling  (hence  the  term  “inner  loop”)  or  using  other  methods  such  as  the  First  Order  Reliability  Method 
(FORM).  Two  primary  cases  where  analytical  methods  can  be  use  used  are  the  normal  format  (linear 
combination  of  normal  random  variables)  and  the  lognormal  format  (multiplicative  combination  of 
lognormal  variables).  The  normal  format  is  used  in  this  paper  to  provide  an  efficient  analysis  method. 

Given  a  response  function  Z  =  A0  +  ^  .4iXi  where,  A0  and  Aj  are  constants  of  the  statistically 

independent  normally  distributed  random  variables  Xr  the  response  mean  and  standard  deviation  can  be 
computed  exactly  as 

Vz=\  +  2  A, .ft  (3) 

and 

<J/  ='  ^jAWi  (4) 

B.  Determination  of  Output  Distribution 

A  schematic  overview  of  the  methodology  is  shown  in  Figure  1.  The  statistical  process  is  explained  as 
follows: 

•  A  distribution  of  each  variable  X  is  provided. 

•  The  population  mean  and  standard  deviation  of  each  random  variable  is  sampled  from  the  t  and 
chi-squared  distribution,  respectively,  see  Eqs.  (1)  and  (2). 

•  The  mean  and  standard  deviation  of  the  response  are  calculated  as  shown  in  Equations  (3)  and 

(4). 

•  Steps  2  and  3  are  repeated  many  times  to  develop  a  population  of  response  moments  fx7  and 

•  Standard  deviation  of  the  response  mean  cr '  is  calculated 

•  An  optimization  algorithm  determines  the  optimal  number  of  additional  experiments  Z)  to 
minimize  the  standard  deviation  of  the  response  moments  subject  to  a  constraint  ^  CiDi  <  B, 
where  C,  is  the  cost  of  each  additional  experiment  and  B  is  the  total  funds  available 

C.  Optimization 

The  constraint  optimization  problem  described  above  is  formulated  as  following: 

Objective:  Minimize  ( a ^  ) 
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Constraint: 


Allocated  resources,  B 


Variable  bounds: 


SN 

,  Cfl,  s  B 


o'  <  D  <  I)' 


There  are  many  optimization  methods  available  to  solve  an  optimization  problem.  The  presence  of 
integer  variables  and  high  computational  cost  associated  with  the  function  evaluation  suggests  that  a 
population-based  approach  is  suitable  for  the  problem.  A  particle  swarm  optimization  (PSO)  modified  to 
handle  integer  variable  is  selected  for  the  problem.  PSO  is  selected  over  other  population-based  approach 
because  of  ease  of  implementation  and  lower  number  of  user  parameters. 

PSO  as  proposed  by  Kennedy  [3],  modified  the  population  from  step-to-step  based  on  a  set  of 
communication  rules.  In  this  set  of  rules,  each  particle  uses  its  current  fitness,  its  best  fitness,  the  best 
fitness  of  all  particles  and  a  communication  structure  to  determine  movement  parameters.  These  parameters 
change  over  the  course  of  iterations,  and  particle  population  tends  to  converge  and  provide  an  optimum 
solution. 

v/+1  =  wvJk  +  <hrji(Pbest j  ~  5?)  +  q2rj2(8best  ~  s ))  (5) 


Equation  (5)  calculates  the  velocity  or  movement  of  each  particle  for  the  k+1  iteration.  This 
equation  uses  velocity  ( v*),  position  ( sk),  personal  best  ( pbest ;  ),  global  best  ( gbest ),  and  three  constants 

(  W ,  qy  q2)  from  iteration  k  to  calculate  the  velocity  of/  particle.  This  velocity  is  used  to  calculate  the  new 
position  of  the  particle  as  shown  in  Equation  (6). 


Sjk+1  =  Skj+v)+1 


(6) 


An  investigation  was  performed  to  determine  if  a  special  strategy  [7]  was  required  to  handle 
integer  variables.  A  number  of  integer  variable  problems  were  solved  to  find  out  if  rounding  the  velocity 
trajectory  of  the  particles  was  sufficient  to  handle  the  integer  variables.  The  results  indicated  that  PSO  was 
able  to  solve  the  problem  and  no  special  strategy  was  required  to  handle  integer  variables.  Manual  search 
was  used  to  determine  appropriate  values  of  the  optimization  variable  for  each  problem  and  once  the 
parameters  were  decided,  these  parameters  were  fixed  for  all  similar  problems.  For  two  and  three  variable 
numerical  problems,  the  PSO  parameters  used  are  given  in  the  following  table. 


Description  (Symbols) 

Value 

Inertia  weight  ( w  ) 

0.529 

Cognitive  constant  ( qt) 

1.494 
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Social  constant  ( q2 ) 


Maximum  velocity 

Du  -Dl 

2 

Minimum  velocity 

Du  -  Dl 

2 

1.494 


The  number  of  particles  and  the  number  of  iterations  were  selected  for  each  problem  depending 
upon  the  number  of  variables  and  complexity  of  the  problem.  For  all  the  problems  in  the  paper,  the 
optimization  used  normalized  variables  (between  0  and  1).  This  normalization  assisted  in  optimizing 
problems  with  variables  of  different  order  of  magnitudes  and  ranges. 

The  formulated  problems  were  constrained  optimization  problems.  The  constraints  were  handled 
using  a  penalty  function  approach.  The  penalty  function  approach  used  here  is  a  modified  form  of  exterior 
penalty  function  method.  In  this  approach,  a  feasible  design  point  objective  function  value  was  always 
better  than  an  infeasible  design  point.  This  was  achieved  by  applying  a  huge  penalty  to  the  objective 
function  value  of  an  infeasible  design  point.  This  strategy  forced  all  the  particles  to  move  toward  the 
feasible  design  space. 

The  optimization  algorithm  employed  two-convergence  criteria  to  stop  the  optimization.  The  first 
criterion  was  that  if  the  best  solution  does  not  change  for  20  consecutive  iterations  then  the  optimization 
stops.  The  second  criterion  was  that  the  optimization  reaches  the  maximum  number  of  iterations  specified 
at  the  beginning  of  the  optimization. 

The  problems  solved  in  the  paper  involve  a  large  number  of  random  numbers  (Monte  Carlo 
samples,  initial  experiments,  and  additional  experiments).  This  randomness  is  also  involved  while 
calculating  the  objective  function  values  during  optimization.  The  randomness  is  a  challenging  issue  to 
obtain  consist  results  from  optimization  because  the  randomness  can  change  the  optimum  solution  each 
time  the  same  problem  is  solved.  This  is  one  of  the  main  conclusions  of  this  study.  As  described  in  the 
numerical  examples,  a  range  of  studies  were  performed  to  investigate  the  effect  of  randomness  in  the 
optimization. 


3.  Numerical  examples 

The  methodology  described  above  was  applied  to  two-  and  three-variable  problems.  All  cases 
considered  a  linear  response  function  and  used  the  normal  format  to  determine  the  response  moments. 

A.  Two-variable  problem 

The  response  was  defined  by  a  linear  function  Z  =  l+10Aj  +  1 0 A,  ,  the  variables  X,  and  X2 
were  normally  distributed  with  /x.  and  cr.and  Ei  original  data  points  were  used.  Each  additional 
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experiment  had  a  cost  of  $1,000  and  the  total  money  available  was  $20,000.  Therefore,  the  maximum 
additional  number  of  experiments  for  each  of  the  variables  was  20. 

A  numerical  study  to  determine  a  sufficient  number  of  Monte  Carlo  samples  was  conducted  using 
case  i  described  below.  The  results  indicated  that  convergence  occurred  at  10,000  samples.  Therefore, 
10,000  samples  were  used  for  all  the  case  studies  described  below. 

In  the  numerical  examples  below,  the  following  notation  is  used. 

•  E,  represents  the  number  of  initial  data  points  for  random  variable  X  , 

•  XE.  represents  the  values  of  the  initial  data  points  for  random  variable  Xjt 

•  Dj  represents  the  number  of  additional  experiments  for  random  variable  X  , 

•  X D  represents  the  values  of  the  additional  experiments  for  random  variable  X r 

•  au  -  orig  represents  the  standard  deviation  of  the  mean  of  Z  based  on  the  original  data  XEj 

•  o '  -  opt  represents  the  optimum  solution  of  the  standard  deviation  of  the  mean  of  Z 

•  XE.  and  the  additional  experiments  X D 
Several  different  cases  were  considered: 

i  E{  =  E2  =10,  XE  *  X E  ,  XD  - random :  That  is,  the  number  of  initial  experiments  were 
set  equal  to  10  ( El  =  E2  =  10),  but  the  initial  values  were  different  ( XE  XE  ).  The  results 
of  the  additional  experiments  ( D:)  were  allocated  randomly  according  to  the  input 
distribution,  e.g.,  normal  with  mean  and  standard  deviation  equal  to  10  and  3,  respectively. 

ii  E{  =  E2  =10,  XE  *  XE  ,  XD  =  XD  =  p'-  Similar  to  case  i  except  the  results  for  any 

additional  experiments  were  set  at  the  mean  values. 

iii  Et  =  E2  =  10,  XE  =  XE  ,  XD  -  random :  Similar  to  case  i  except  the  initial  experimental 
results  were  the  same  and  allocated  randomly  according  to  the  input  distributions 

iv  El  =  E,  =10,  XE  =  XE  ,  XD  =  X D  =  U :  The  values  for  the  initial  experiments  are  the 
same  and  the  results  of  any  new  experiments  set  at  the  mean  of  the  distributions. 

v  Et  =  20,  E2  =  10,  XE  *  X ,  XD  -  random :  A  different  number  of  initial  data  points  and 

values  were  used  for  Xt  and  X2-  The  additional  experiments  were  allocated  randomly 
according  to  the  input  distribution. 

vi  Et  =20,  E2  =  10,  XE  ^  XE  ,  XD  =  XD  =  p'.  Similar  to  case  v  except  the  values  for  the 
additional  experiments  were  located  at  mean  of  the  variable. 

Case  i:  Ex  =  E2  =  10,  XE  ^  XE  ,  XD  -  random 

The  parameters  are  defined  in  Table  1  and  the  design  space  is  shown  in  Figure  2;  the  colored  area 
indicates  the  feasible  area  while  the  flat  blue  area  represents  the  infeasible  area.  The  values  for  the  initial 
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data  points,  XE  ,  are  shown  in  Table  2.  The  values  for  the  results  of  the  additional  experiments,  X n  ,  were 
determined  randomly  according  to  the  random  variable  distribution.  For  example,  the  values  chosen  for 
XD  were  obtained  randomly  from  a  normal  distribution  with  mean  and  standard  deviation  of  10  and  3, 
respectively.  Because  of  this  randomness,  different  results  were  obtained  for  Z)  for  each  analysis. 

The  optimization  process  was  conducted  using  10  iterations  and  20  particles.  Each  addition 
experiment  cost  $1000  and  the  total  amount  funds  available  was  $20,000.  The  particles  during  the  1st,  5th 
and  10th  iteration  for  one  particular  analysis  is  shown  in  Figure  3.  Table  3  show  the  results  obtained  from 
running  the  analysis  5  separate  times.  The  effect  of  the  randomness  in  the  values  of  XE  and  XD  is  evident; 
a  large  variation  in  the  optimal  number  of  additional  experiments  for  each  variable  was  obtained.  The 
problem  appears  symmetric  with  respect  to  Xx  and  X 2',  therefore,  the  intuitive  solution  is  Dl  =  D,  =10. 
However,  due  to  the  random  nature  of  the  actual  values  used  during  the  analysis,  quite  different  answers 
can  be  obtained. 

In  all  cases,  the  additional  experiments  reduced  the  standard  deviation  of  jl  by  approximately 

30%.  For  example,  in  run  1,  the  standard  deviation  of  the  response  mean  without  any  additional  data  was 
18.28  and  was  reduced  to  11.23  after  adding  the  additional  experiments. 

The  PDF  of  the  response  mean  is  shown  in  Figure  4  for  run  1.  The  red  dotted  line  shows  the  PDF  of  the 
response  mean  previous  to  any  additional  data,  the  solid  blue  line  shows  the  PDF  after  adding  additional 
experiments.  A  reduction  in  the  standard  deviation  of  jl_  is  evident.  Note,  also  that  due  to  the  random 

allocation  of  the  results  from  the  additional  experiments,  the  mean  after  adding  the  additional  experiments 
can  shift  from  the  original  data. 


Table  1.  Variable  Parameters  (Case  i) 


A- 

/Zy 

<Jx, 

E, 

C, 

*1 

10 

10 

3 

10 

$1,000 

x2 

10 

10 

3 

10 

$1,000 

Table  2.  Values  for  variables  XE  and  XD  (Case  i) 


9.14 

11.71 

9.18 

11.57 

9.04 

10.76 

11.19 

8.83 

8.61 

11.86 

10.52 

7.85 

2.05 

8.98 

8.46 

6.98 

5.39 

6.94 

6.88 

5.81 

9.31 

9.05 

5.13 

8.95 

xEl 

9.81 

9.25 

14.33 

9.41 

5.89 

11.10 

5.38 

14.44 

7.09 

7.26 

8.60 

10.04 

10.78 

8.73 

11.80 

9.34 

16.06 

11.31 

12.62 

11.47 
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5.91 

12.35 

7.62 

12.30 

10.86 

11.16 

11.32 

11.41 

7.86 

8.44 

8.20 

11.43 

11.55 

11.04 

5.84 

11.43 

Table  3.  Summary  (Case  i) 


Run 

A 

d2 

-orig 

°,7  -opt 

%  Reduction 

2C,Z> 

1 

15 

5 

18.28 

11.23 

38.6% 

$20,000 

2 

16 

4 

20.46 

13.68 

33.1% 

$20,000 

3 

6 

14 

19.18 

12.69 

33.8% 

$20,000 

4 

11 

9 

21.61 

14.83 

31.4% 

$20,000 

5 

7 

13 

25.17 

16.76 

33.4% 

$20,000 

Case  ii:  E,  =  E2  =  10,  XEi  *  X ^  XD]  =  X,h  =  n 

For  this  case,  the  additional  experiments,  XD  were  set  the  mean  of  each  variable,  that  is,  it  was 

assumed  that  after  conducting  the  new  experiments,  the  results  were  equal  to  the  population  mean.  Clearly, 
this  is  an  assumption  not  realized  in  practice  but  useful  to  examine  the  affect  of  randomness  in  the 
optimization  process. 

The  design  space  is  shown  in  Figure  5,  the  colored  area  shows  the  feasible  area,  while  the  flat 
blue  area  represents  the  infeasible  area.  The  surface  shown  in  Figure  5  is  smoother  than  in  case  i  due  to  the 
removal  of  some  randomness.  The  optimization  process  used  10  iterations  with  10  particles.  PSO  require 
lower  number  of  particles  as  compared  to  Case  i  because  the  removal  of  randomness  provides  a 
monotonous  (smooth)  design  space. 

The  results  for  5  analyses  are  shown  in  Table  6.  The  results  are  a  little  more  symmetric  than  case  i,  as 
expected,  but  still  show  a  significant  amount  of  asymmetry.  The  optimization  converged  after  a 
combination  of  10  iterations  and  10  particles  and  the  optimization  process  is  represented  in  Figure  6.  The 
addition  of  the  new  experimental  data  reduced  the  standard  deviation  of  the  response  mean  by 
approximately  50%  in  each  run.  The  PDF  of  the  response  mean  is  shown  in  Figure  7.  The  red  dotted  line 
shows  the  PDF  of  the  response  mean  previous  to  any  additional  data,  the  solid  blue  line  shows  the  PDF 
after  adding  a  new  set  of  data. 


Table  4.  Variable  Parameters  (Case  ii) 


A- 

Ay, 

A 

C, 

A 

10 

10 

3 

10 

$1,000 

x2 

10 

10 

3 

10 

$1,000 

Table  5.  Values  for  variables  XE  and  XD  (Case  ii) 


10 


*D2 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

Xe2 

9.81 

9.25 

14.33 

9.41 

5.89 

11.10 

5.38 

14.44 

7.09 

7.26 

8.60 

10.04 

10.78 

8.73 

11.80 

9.34 

16.06 

11.31 

12.62 

11.47 

Table  6.  Summary  (Case  ii) 


Run 

A 

d2 

-orig 

°u7  ~°P> 

%  Reduction 

2C,D, 

1 

13 

1 

18.28 

8.51 

53.4% 

$20,000 

2 

13 

7 

20.46 

10.43 

49.0% 

$20,000 

3 

9 

11 

19.18 

9.37 

51.1% 

$20,000 

4 

10 

10 

21.61 

11.33 

47.6% 

$20,000 

5 

11 

9 

25.17 

12.61 

49.9% 

$20,000 

Case  iii:  Et  =  E2  =  10,  XE  =  XE,  XD  -  random 

For  this  case,  the  initial  values  of  the  variables  are  equal,  XE  =  XE  ,  and  the  values  of  the  additional 
experiments  were  determined  randomly  according  to  the  random  variable  distributions.  The  problem 
parameters  ( An  jlx  ,  ox  ,  Er  C()  are  shown  in  Table  7  and  the  experimental  location  values  ( XE  ,  XD  ) 

are  shown  in  Table  8.  The  design  space  is  shown  in  Figure  8  where  the  colored  area  shows  the  feasible 
area,  while  the  flat  blue  area  represents  the  infeasible  area.  10  iterations  and  20  particles  were  used  during 
the  optimization  process.  The  particle  movement  for  first,  fifth  and  tenth  iteration  is  shown  in  Figure  9. 

The  asymmetry  in  the  allocation  of  the  additional  experiments  seemed  to  be  reduced  compared  to  both 
cases  i  and  ii;  however,  more  analyses  are  needed  to  confirm  this  conclusion.  The  addition  of  the  new 
experimental  data  reduced  the  standard  deviation  of  the  response  mean  by  close  to  40%. 

The  PDF  of  the  response  mean  is  shown  in  Figure  10.  The  red  dotted  line  shows  the  PDF  of  the 
response  mean  previous  to  any  additional  data,  the  solid  blue  line  shows  the  PDF  after  adding  a  new  set  of 
data. 
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Table  7.  Variable  Parameters  (Case  iii) 


A- 

^x, 

<*x, 

E, 

C, 

*1 

10 

10 

3 

10 

$1,000 

A 

10 

10 

3 

10 

$1,000 

Table  8.  Values  for  variables  XE  and  XD  (Case  iii) 


A, 

A2 

12.51 

8.04 

7.83 

6.76 

7.84 

9.86 

9.40 

11.14 

9.94 

9.01 

10.84 

8.50 

13.17 

9.89 

11.87 

9.48 

4.75 

7.13 

12.09 

13.88 

12.43 

11.32 

11.91 

13.84 

13.93 

8.51 

10.98 

6.64 

7.98 

12.42 

9.55 

10.12 

2.65 

7.73 

11.42 

9.73 

10.35 

3.97 

8.23 

13.25 

At, 

a2 

9.96 

9.96 

8.26 

8.26 

16.41 

16.41 

9.23 

9.23 

5.77 

5.77 

15.31 

15.31 

10.98 

10.98 

6.64 

6.64 

11.86 

11.86 

13.81 

13.81 

Table  9.  Summary  (Case  iii) 


Run 

A 

A 

-orig 

°,7  ~°P< 

%  Reduction 

2C,.Z> 

1 

8 

12 

23.34 

13.53 

42.0% 

$20,000 

2 

14 

6 

19.61 

13.53 

31.0% 

$20,000 

3 

10 

10 

22.88 

14.26 

37.7% 

$20,000 

4 

9 

11 

21.47 

13.12 

38.9% 

$20,000 

5 

10 

10 

22.88 

14.26 

37.7% 

$20,000 

Case  iv:  E}  =  E2  =  10,  XEi  =  XEi,  XD  =  XD[  =  [i 

For  this  case,  the  initial  values  of  the  variables  were  equal,  XE  =  XE  ,  and  the  values  of  the  additional 

experiments  were  also  equal  and  set  at  the  mean  value,  XD  =  XD  =  Ll  determined  according  to  the 

random  variable  distribution.  The  parameters  are  shown  in  Table  10  and  Table  11.  10  iterations  and  10 
particles  were  used.  The  design  space  is  shown  in  Figure  1 1  where  the  colored  area  shows  the  feasible  area, 
while  the  flat  blue  area  represents  the  infeasible  area. 
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The  additional  experiments  of  each  variable  were  10  and  10  respectively,  as  expected;  in  this  case,  there 
is  no  asymmetry.  The  optimization  converged  after  a  combination  of  10  iterations  and  10  particles  and  the 
optimization  process  is  represented  in  Figure  12.  The  addition  of  the  new  experimental  data  reduced  the 
standard  deviation  of  the  response  mean  from  a  value  of  23.34  to  a  value  of  11.53;  a  reduction  of  50.6%. 

The  PDF  of  the  response  mean  is  shown  in  Figure  13.  The  red  dotted  line  shows  the  PDF  of  the 
response  mean  previous  to  any  additional  data,  the  solid  blue  line  shows  the  PDF  after  adding  a  new  set  of 
data.  The  standard  deviation  of  the  response  mean  was  reduced  by  approximately  50%. 


Table  10. _ Variable  Parameters  (Case  iv) 


A- 

Px, 

E, 

C, 

*1 

10 

10 

3 

10 

$1,000 

x2 

10 

10 

3 

10 

$1,000 

Table  11.  Values  for  variables  XE  and  XD  (Case  iv) 


10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

9.96 

9.96 

8.26 

8.26 

16.41 

16.41 

9.23 

9.23 

5.77 

5.77 

15.31 

15.31 

10.98 

10.98 

6.64 

6.64 

11.86 

11.86 

13.81 

13.81 

Table  12. _ Summary  (Case  iv) 


Run 

A 

D2 

<Tt,z  -orig 

°,z-°Pt 

%  Reduction 

2C,D, 

1 

10 

10 

23.34 

11.53 

50.6% 

$20,000 

2 

10 

10 

14.86 

7.3 

50.9% 

$20,000 

3 

10 

10 

15.38 

7.87 

48.8% 

$20,000 

4 

10 

10 

12.36 

6.05 

51.1% 

$20,000 

5 

10 

10 

16.96 

8.41 

50.4% 

$20,000 
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Ccisev:  E{=  20,  E2=  10,  X  E  *  X  E,  X n  -random 

For  this  case,  the  number  of  initial  experiments  and  the  values  of  the  experiments  were  not  equal;  and 
the  numbers  of  initial  data  points  were  20  and  10  respectively,  and  the  values  of  the  additional  experiments 
were  determined  according  to  the  random  variable  distributions.  The  parameters  are  shown  in  Table  13  and 
the  design  space  is  show  in  Figure  14  where  the  colored  area  shows  the  feasible  area,  while  the  flat  blue 
area  represents  the  infeasible  area.  The  optimization  used  10  iterations  and  20  particles.  The  optimization 
process  is  represented  in  Figure  15  for  run  1. 

In  this  case,  the  number  of  additional  experiments  for  each  variable  was  heavily  weighted  towards  more 
experiments  of  X2,  see  Table  15.  This  occurs  largely  because  more  information  is  already  available  for 
Xp  therefore,  the  variance  of  the  population  mean  jdl  is  tighter  than  jd2-  The  addition  of  the  new 
experimental  data  reduced  the  standard  deviation  of  the  response  mean  by  20-40%. 

The  PDF  of  the  response  mean  is  shown  in  Figure  16.  The  red  dotted  line  shows  the  PDF  of  the 
response  mean  previous  to  any  additional  data,  the  solid  blue  line  shows  the  PDF  after  adding  a  new  set  of 
data. 


Table  13.  Variable  Parameters  (Case  v) 


A 

/% 

°xx 

E, 

c, 

*1 

10 

10 

3 

20 

$1,000 

X2 

10 

10 

3 

10 

$1,000 

Table  14.  Values  for  variables  XE  and  XD  (Case  v) 


Xe2 

11.60 

7.91 

14.10 

18.44 

7.73 

12.86 

12.34 

7.45 

9.31 

10.81 

12.46 

12.02 

7.68 

12.16 

8.00 

10.75 

13.68 

9.31 

8.54 

13.77 

4.72 

12.68 

9.32 

6.81 

10.23 

9.01 

4.53 

5.86 

12.73 

7.61 

11.80 

8.03 

5.29 

8.32 

13.69 

12.63 

9.20 

8.00 

6.38 

10.83 

7.17 

4.11 

13.95 

7.40 

5.03 

9.37 

9.63 

5.62 

9.03 

5.13 

16.91 

9.96 

9.32 

14.06 

12.42 

9.70 

9.30 

15.27 

9.25 

12.82 

12.60 

9.31 

13.45 

10.98 

4.78 

7.25 

9.70 

11.78 

8.78 

10.11 
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Table  15. _ Summary  (Case  v) 


Run 

A 

d2 

-orig 

°,7  -°Pl 

%  Reduction 

2C,D, 

1 

l 

19 

16.99 

12.49 

26.5% 

$20,000 

2 

3 

17 

13.96 

11.18 

19.9% 

$20,000 

3 

1 

18 

15.77 

10.38 

34.2% 

$19,000 

4 

6 

13 

17.94 

12.05 

32.8% 

$19,000 

5 

1 

19 

19.98 

11.54 

42.2% 

$20,000 

Casevi:  Et  =20,  E2  =10,  XEi  *  XE  ,  XD  =  XD[  =  /U 

For  this  case,  the  number  of  initial  experiments  and  the  values  of  the  experiments  were  not  equal,  the 
numbers  of  initial  data  points  were  20  and  10  respectively,  and  the  values  of  the  additional  experiments 
were  set  at  the  mean  value  of  the  random  variables.  The  parameters  are  shown  in  Table  16  and  the  design 
space  is  show  in  Figure  17  where  the  colored  area  shows  the  feasible  area,  while  the  flat  blue  area 
represents  the  infeasible  area.  For  this  case,  the  number  of  additional  experiments  were  weighted  towards 
X2  but  not  as  heavily  as  in  case  v,  see  Table  18.  The  optimization  converged  after  a  combination  of  10 
iterations  and  10  particles  and  the  optimization  process  is  represented  in  Figure  18.  The  addition  of  the  new 
experimental  data  reduced  the  standard  deviation  of  the  response  mean  by  approximately  40%. 

The  PDF  of  the  response  mean  is  shown  in  Figure  19.  The  red  dotted  line  shows  the  PDF  of  the 
response  mean  previous  to  any  additional  data,  the  solid  blue  line  shows  the  PDF  after  adding  a  new  set  of 
data. 


Table  16.  Variable  Parameters  (Case  vi) 


A- 

A%( 

<*x, 

E, 

C, 

10 

10 

3 

20 

$1,000 

*2 

10 

10 

3 

10 

$1,000 

Table  17.  Values  for  variables  XE  and  XD  (Casevi) 


** 

11.60 

7.91 

14.10 

18.44 

7.73 

12.86 

12.34 

7.45 

9.31 

10.81 

12.46 

12.02 

7.68 

12.16 

8.00 

10.75 

13.68 

9.31 

8.54 

13.77 

4.72 

12.68 

9.32 

6.81 

10.23 

9.01 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

15 


4.53 

5.86 

12.73 

7.61 


10 

10 

10 

10 

10 

10 

10 

10 

Table  18. _ Summary  (Case  vi) 


Run 

D\ 

d2 

°rz  ~oriS 

-°Pf 

%  Reduction 

2C,D, 

1 

6 

14 

16.99 

9.6 

43.5% 

$20,000 

2 

13 

7 

13.34 

8.76 

34.3% 

$20,000 

3 

6 

14 

19.28 

10.43 

45.9% 

$20,000 

4 

8 

12 

21.1 

11.93 

43.5% 

$20,000 

5 

7 

13 

18.49 

10.16 

45.1% 

$20,000 

B.  Three  random  variables 

The  response  was  define  by  a  linear  function  Z  =  1+  7Z,  +  3Z^  +  13Z,,  the  variables  Xr  X2  and  X3 
were  normally  distributed.  Different  values  were  used  for  A’s,  input  PDF  parameters,  initial  data  points  and 
cost  functions.  The  total  funds  available  was  $20,000.  Two  different  cases  were  considered: 

i  The  results  of  the  additional  experiments  were  allocated  randomly  according  to  the 
input  distribution. 

ii  The  results  of  the  additional  experiments  are  allocated  at  the  mean  of  the  variable. 

Case  vii:  £,  =  5,  E,  =  1 5,  £,=10,  XE^  *XE  *XE^,  XD  -random 

The  parameters  are  described  in  Table  19  and  Table  20,  the  values  for  the  results  of  the  additional 
experiments  were  determined  according  to  the  random  variable  distribution. 


Table  19.  Variable  Parameters  (Case  vii-Three  variables) 


A 

<*x, 

E, 

c, 

i 

10 

3 

5 

$700 

X2 

3 

0 

1 

15 

$500 

x3 

13 

7 

2 

10 

$1,000 

Table  20.  Values  for  variables  XE  and  XD  (Case  vii-Three  variables) 


9.07 

-0.35 

9.07 

11.64 

-0.58 

6.57 

13.93 

0.14 

5.12 

14.92 

-0.79 

7.10 

14.69 

-0.14 

10.62 

-0.15 

7.15 

^3 

9.30 

1.51 

9.30 

7.33 

-0.52 

9.68 

11.75 

1.75 

10.32 

12.61 

1.10 

9.79 

14.71 

-0.02 

8.51 

8.47 

-0.28 

6.67 
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-0.21 

7.18 

-0.41 

5.68 

0.15 

9.15 

0.47 

7.51 

1.76 

-0.41 

1.00 

0.56 

-0.57 

12.40 

1.33 

6.18 

11.85 

0.46 

4.16 

12.21 

-1.42 

6.84 

8.89 

0.39 

6.57 

12.69 

0.15 

7.51 

9.53 

1.06 

7.46 

9.05 

1.03 

7.52 

8.77 

-0.40 

6.19 

11.43 

-1.16 

6.64 

15.48 

-1.80 

6.28 

8.97 

-0.17 

5.99 

9.44 

0.88 

8.98 

4.78 

0.29 

9.98 

11.26 

0.76 

5.16 

11.93 

-2.03 

12.80 

0.04 

8.18 

-0.13 

8.97 

0.40 

7.45 

1.47 

13.13 

-0.10 

11.29 

-0.73 

10.40 

1.52 

-0.45 

-2.05 

-1.04 

-1.74 

0.62 

0.47 

-0.24 

0.31 

-0.52 

-0.34 

-0.26 

1.43 

The  following  study  was  performed  to  investigate  the  behavior  of  the  problem.  The  numbers  of  additional 
experiments  were  determined  using  a  simple  allocation  scheme.  For  example,  equal  money  for  each 
variable  or  maximum  allocation  for  each  variable,  the  summary  is  shown  in  Table  21.  Thus,  Table  21  show 
the  results  that  can  be  expected  without  optimization. 


Table  21. _ Summary  of  study  (Case  vii-Three  variables) 


D, 

28 

0 

0 

9 

9 

d2 

0 

40 

0 

13 

9 

O, 

0 

0 

20 

7 

9 

°U2  -°Pf 

10.87 

15.34 

12.56 

10.59 

10.77 

%  Reduction 

30% 

1% 

19% 

31% 

30% 

2C,A 

19600 

20000 

20000 

19800 

19800 
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The  results  from  the  optimization  process  are  shown  in  Table  22.  The  standard  deviation  of  the  response 
mean  was  reduced  by  typically  40%.  The  PDF  of  the  response  mean  is  shown  in  Figure  20  for  run  1.  20 
iterations  and  40  particles  were  used.  The  red  dotted  line  shows  the  PDF  of  the  response  mean  previous  to 
any  additional  data,  the  solid  blue  line  shows  the  PDF  after  adding  a  new  set  of  data.  Due  to  the  random 
allocation  of  the  results  from  the  additional  experiments,  the  mean  after  adding  the  additional  experiments 
shifted  to  the  right. 


Table  22,  Summary  (Case  vii-Three  variables) 


Run 

A 

A 

A 

An  -orig 

°,z-°Pt 

%  Reduction 

2C,D, 

1 

18 

0 

7 

15.45 

9.61 

38% 

$19,600 

2 

10 

0 

13 

19.25 

10.06 

48% 

$20,000 

3 

11 

0 

12 

26.22 

14.03 

46% 

$19,700 

4 

17 

0 

8 

21.62 

11.12 

49% 

$19,900 

5 

9 

1 

13 

13.67 

11.39 

17% 

$19,800 

Whit  the  same  parameters  set  in  Table  22,  the  standard  deviation  of  the  response  mean  is  calculated  by 
changing  the  additional  number  of  experiments  by  a  small  number  and  compared  to  the  optimum  solution 
of  9.61.  This  study  was  performed  to  determined  the  robustness  of  the  solution.  Table  23  shows  different 
possibilities. 


Table  23. _  Study  around  the  optimal  solution  (Case  vii-Three  variables) 


D, 

17 

18 

17 

16 

15 

d2 

0 

0 

0 

0 

0 

O, 

7 

6 

8 

8 

9 

°,z-°pt 

9.72 

9.74 

9.86 

10.22 

9.84 

%  Increment 

1% 

1% 

3% 

6% 

2% 

2C,D, 

18900 

18600 

19900 

19200 

19500 

Case  viii:  Et  =5,  E2=  15,  E 3  =  10,  X£|  ^  XE ^  *  XE};  XDi  =  10,  XD  =  0 ,  XD  =  7 

The  parameters  are  described  in  Table  24,  the  values  for  the  results  of  the  additional  experiments  were 
set  at  the  mean.  The  optimization  process  obtained  an  optimal  number  of  additional  experiments  of  as 
shown  in  Table  26.  20  iterations  and  40  particles  were  used.  The  standard  deviation  of  the  response  mean 
was  reduced  by  approximately  60%.  The  PDF  of  the  response  mean  is  shown  in  Figure  21.  The  red  dotted 
line  shows  the  PDF  of  the  response  mean  previous  to  any  additional  data,  the  solid  blue  line  shows  the  PDF 
after  adding  a  new  set  of  data 


Table  24.  Variable  Parameters  (Case  viii-Three  variables) 


A, 

Ay, 

A 

c, 

A 

7 

10 

3 

5 

$700 

*2 

3 

0 

1 

15 

$500 

A 

13 

7 

2 

10 

$1,000 
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4.  Conclusions 


Several  simple  probabilistic  models  were  explored  to  assess  the  optimum  additional  experiments  to 
reduce  the  output  variance  in  the  mean  or  standard  deviation  of  a  probabilistic  response.  The  example 
problems  were  linear  combinations  of  normal  random  variables  such  that  the  inner  loop  could  be  computed 
in  closed  form  to  facilitate  efficient  exploration;  however,  methodology  is  not  limited  to  such  simple 
problems.  An  efficient  optimization  strategy  was  implemented  using  a  particle  swarm  optimization  (PSO) 
method  modified  to  handle  integer  variables.  PSO  was  selected  over  other  population-based  approaches 
because  of  ease  of  implementation  and  a  lower  number  of  user  parameters. 

The  number  of  additional  experiments  to  add  for  each  random  variable  to  reduce  the  standard  deviation 
of  the  response  mean  depend  upon  many  factors:  the  number  of  initial  data  points  (  E  ),  the  importance  of 
the  random  variable  in  the  response  (A,),  the  range  of  the  mean  of  the  random  variable  ( jlt),  the  cost  of 
each  experiment  ( C.).  In  addition,  the  actual  values  affect  the  solution.  For  complicated  problems,  it  is  not 
possible  to  ascertain  an  optimal  allocation  of  additional  experiments,  hence  the  value  of  the  approach 
described  here. 

The  results  indicate  that  the  optimization  is  difficult  due  primarily  to  the  randomness  in  the  results 
of  the  additional  experiments.  The  optimization  results  were  dependent  upon  how  the  additional 
experiments  were  chosen.  If  chosen  randomly,  the  optimum  experiments  were  also  random,  that  is, 
different  results  were  obtained  if  the  same  analysis  was  rerun;  however,  general  trends  were  observed. 
More  consistent  results  were  obtained  if  the  additional  experiments  were  all  located  at  the  mean  value  of 
the  input  random  variables.  This  was  evident  in  the  symmetric  two  random  variable  problem  where  the 
intuitive  results  that  and  equal  number  of  both  experiments  were  needed  would  often  not  be  obtained  unless 
the  initial  and  additional  experimental  results  were  selected  to  be  the  same  for  both  variables.  An 
examination  of  the  three  random  variable  problem  clearly  showed  the  value  of  the  optimization  algorithm 
since  a  larger  reduction  in  the  standard  deviation  of  the  response  mean  was  obtained  compared  to  ad  hoc 
approaches. 
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6.  Figures 


Figure  1.  Schematic  flow  chart  of  the  optimal  allocation  methodology 
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Figure  2.  Design  Space  (Two  variable  -  Case  i) 


Figure  3.  Particle  movement  (Two  variable  -  Case  i)  a)  Particles  at  first  iteration,  b)  Particle  at 
fifth  iteration  and  c)  Particles  at  the  last  iteration  where  the  blue  dot  represents  the  optimum  solution. 
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Figure  4.  PDF  of  Response  mean  (Two  variables  -  Case  i) 
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Figure  5.  Design  Space  (Two  variable  -  Case  ii) 
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Figure  6.  Particle  movement  (Case  ii)  a)  Particles  at  first  iteration,  b)  Particle  at  fifth  iteration  and 
c)  Particles  at  the  last  iteration  where  the  blue  dot  represents  the  optimum  solution. 


Figure  7.  PDF  of  Response  mean  (Two  variables  -  Case  ii) 
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Figure  8.  Design  Space  (Two  variable  -  Case  iii) 


Figure  9.  Particle  movement  (Case  iii)  a)  Particles  at  first  iteration,  b)  Particle  at  fifth  iteration 
and  c)  Particles  at  the  last  iteration  where  the  blue  dot  represents  the  optimum  solution. 
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Figure  10.  PDF  of  Response  mean  (Two  variables  -  Case  iii) 
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Figure  11.  Design  Space  (Two  variable  -  Case  iv) 


28 


PDF 


Figure  12.  Particle  movement  (Case  iv)  a)  Particles  at  first  iteration,  b)  Particle  at  fifth  iteration 
and  c)  Particles  at  the  last  iteration  where  the  blue  dot  represents  the  optimum  solution. 


Figure  13.  PDF  of  Response  mean  (Two  variables  -  Case  iv) 
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Figure  14. 


Design  Space  (Two  variable  -  Case  v) 


Figure  15.  Particle  movement  (Case  v)  a)  Particles  at  first  iteration,  b)  Particle  at  fifth  iteration  and 
c)  Particles  at  the  last  iteration  where  the  blue  dot  represents  the  optimum  solution. 


30 


Standard  Deviation  PDF 


Figure  16.  PDF  of  Response  mean  (Two  variables  -  Case  v) 
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Figure  17.  Design  Space  (Two  variable  -  Case  vi) 
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Figure  18.  Particle  movement  (Case  vi)  a)  Particles  at  first  iteration,  b)  Particle  at  fifth  iteration 
and  c)  Particles  at  the  last  iteration  where  the  blue  dot  represents  the  optimum  solution. 


Figure  19.  PDF  of  Response  mean  (Two  variables  -  Case  vi) 
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Figure  20.  PDF  of  the  Response  mean  (Three  variables  -  Case  vii) 


Figure  21.  PDF  of  the  Response  mean  (Three  variables  -  Case  viii) 
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